o 



X 



Theoretical analysis of Sine-collocation methods and Sinc-Nystrom 

methods for initial value problems 

Tomoaki Okayama 

Graduate School of Economics, Hitotsubashi University, 2-1 Naka, Kunitachi, Tokyo 186-8601, Japan 



Abstract 

A Sinc-coUocation method has been proposed by Stenger, and he also gave theoretical analysis of 



u 

^ ■ the method in the case of a 'scalar' equation. This paper extends the theoretical results to the case 

. of a 'system' of equations. Furthermore, this paper proposes more efficient method by replacing 

^ ! the variable transformation employed in Stenger's method. The efficiency is confirmed by both 
of theoretical analysis and numerical experiments. In addition to the existing and newly-proposed 

<;^ ■ Sinc-coUocation methods, this paper also gives similar theoretical results for Sinc-Nystrom meth- 

^ ■ ods proposed by Nurmuhammad et al. From a viewpoint of the computational cost, it turns out 

_H ■ ffiat the newly-proposed Sinc-coUocation method is the most efficient among those methods. 

"^ • 
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i/-^ ' 1. Introduction 

^' 

T^ ' The concem of this paper is a system of initial value problems of the form 

m ■ \y'it) = K{t)y{t) + g{t), a<t<b. 




r. 



(1) 




where K(t) is an n x n matrix, and y(t), g{t), r are column vectors of order n. For Eq. 
eral numerical methods based on the Sine approximation have been developed so far [|2|,J|5, 
^.' and in general, those methods converge exponentially. For example, Carlson et al. PO pro- 
posed a Sinc-coUocation method for Eq. dT), and they also claimed that its convergence rate is 
0{N^ exp(-c VA^)). However, their method is designed not for the finite interval [a, b\ but for the 
infinite interval (-oo, oo) or [0, £»), and users have to know solution's behavior as r ^ oo to imple- 
ment the method. In addition, users also have to know solution's regularity for implementation. 
It is not so practical to assume that solution's behavior and regularity can be known in prior to 
computation, since the solution is an unknown function. 
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Instead of solving Eq. ([B, Stenger [Il3h firstly transformed the problem to the Volterra integral 
equation of the second kind: 



y(t) = r+ {g{s) + K(s)y(s)} ds, a < t < b, (2) 

Ja 

and derived a Sinc-coUocation method for Eq. (|2]). His method does not require solution's be- 
havior as f ^ oo when the given interval [a, b] is finite. In addition, he showed theoretically that 
solution's regularity needed for implementation can be found from the known functions. This is an 
advantage of his method over that of Carlson et al. Moreover, he also showed that the convergence 
rate of his method is 0( VA^exp(-c VA^)), where c is the same constant as the result of Carlson et 
al. It should be noted that those theoretical results were shown only in the case where Eq. <^ is 
a scalar equation (n = 1), although the method was proposed for a system of equations. This is 
because the analysis relies on the explicit form of the solution y that holds only in the scalar case. 

The first objective of this study is to extend Stenger's theoretical results to a system of equa- 
tions. That is, this paper shows that even in the case of a system of equations, solution's regularity 
actually can be found from the known functions K{t) and g{t), and also shows that his method 
converges with the rate: 0( VA^exp(-c ViV)). 

The second objective, which is more important in this paper, is to improve Stenger's method. 
The main idea here is replacement of the variable transformation; the "Single-Exponential trans- 
formation" (SE transformation) is employed in the method of Stenger (and also Carlson et al.), 
but it is replaced with the "Double-Exponential transformation" (DE transformation) in the pro- 
posed method. Those two methods are referred to as the SE-Sinc-coUocation method and the 
DE-Sinc-coUocation method, respectively. It has been known that the replacement of the variable 
transformation often accelerates the convergence flSl [l5||, and in fact, this paper shows by theoret- 



ical analysis that the rate is drastically improved to 0(exp(-cW/ log A^^)) by the replacement. 

The third objective of this study is to give similar theoretical results for Sinc-Nystrom methods 
for Eq. (|2l). The methods have been proposed by Nurmuhammad et al. [5], where both the SE 
transformation and the DE transformation are considered. Those two methods are referred to as 
the SE-Sinc-Nystrom method and the DE-Sinc-Nystrom method, respectively. Any convergence 
analysis has not been given for both Sinc-Nystrom methods, and as it stands users have no clue 
to decide which to choose out of the four methods: SE/DE-Sinc-coUocation methods and SE/DE- 
Sinc-Ny Strom methods. To improve the situation, this paper analyzes the errors theoretically, and 
shows that the convergence rate of the SE-Sinc-Nystrom method is 0(exp(-c VA^)), and that of 
the DE-Sinc-Nystrom method is Oi^f exp{-c'N/ logN)). 

From a viewpoint of the convergence rate, the DE-Sinc-Nystrom method seems to be the best 
among the four methods. From a viewpoint of the computational cost, however, the DE-Sinc- 
coUocation method has several advantages (see discussion in Section 1431) . Moreover, according 
to the theoretical analysis in this study, the difference of the convergence rate between the two is 
quite small, and in fact we can confirm it in numerical experiments (see Section[5]). Therefore, the 
(proposed) DE-Sinc-coUocation method compares favorably with the DE-Sinc-Nystrom method. 

The remainder of this paper is organized as follows. In Section |21 basic definitions and theo- 
rems of Sine methods are stated. In Section |3l four numerical methods to be considered: SE/DE- 
Sinc-Nystrom methods and SE/DE-Sinc-coUocation methods are described. Main theoretical re- 



suits are stated in Section |4l and their proofs are given in Section [6l Numerical examples are 
presented in Section |5l 

2. Basic definitions and theorems of Sine methods 

In this section, fundamental approximation formulas derived from the Sine approximation are 
explained with their convergence theorems. 

2.1. Sine approximation and Sine indefinite integration over the real axis 
The Sine approximation is expressed as 

N 

F(x) « Yj ^0"/^)5 U, h)(x), xeR, (3) 

j=-N 

where the basis function 5(j, h){x) (the so-called Sine function) is defined by 

smnjx/h-j) 

Sij,h)ix) = — — — , 

K{x/h - ]) 

and /z is a step size appropriately chosen depending on a given positive integer N. The Sine 
indefinite integration is derived by integrating both sides of Eq. © as 

F(t) dt^Yj ^0'/^) S U, h){t) dt=Yj FUh)Ja h)ix), x e M, (4) 

-'-~ j=^N ^-°° j=-N 

where 7(j, h){x) is defined by 

7(7, h)(x) = h\^ + - Si[K(x/h - j)] 

Here, Si(A:) is the so-called "sine integral" function, whose routine is available in some numerical 
libraries (IMSL, NAG, GSL, and so on). The approximation dD is called the Sine indefinite 
integration. 

2.2. (Generalized) SE-Sinc approximation and SE-Sinc indefinite integration 

When the target interval [a, b] is finite, the Single-Exponential (SE) transformation 

b - a . I x\ b + a 



c^p b - a , /x\ 

r = r(.) = ^tanh(-)+ 2 ' 

x = r(0 = {^^V(0 = log(^) 



is frequently used with the formulas <^ and (H). Since this transformation maps t 6 (a, b) onto 

x 6 M, we can use Q and © as 

N 

fit) = if o r")icf>'\t)) « Yj f(^P^ o; h)i<f>'^t)), (5) 

j=-N 

fis)ds= /(^^^^(x)){^^^}'(x)dx« 2]/(^f){'A"}'o-/^)^(7;/^)(0"(O), (6) 

Ja J-oo j=-A' 



where t^ = iff^^(jh). The approximations © and © are called the SE-Sinc approximation and the 
SE-Sinc indefinite integration, respectively. 

If / is non-zero at the endpoints t = a and t = b, the SE-Sinc approximation does not work 
accurately near the endpoints, because the right hand side of © tends to when t ^> a and t ^> b. 
To remedy the issue, Stenger [Il3l] introduced the auxiliary basis functions wdt) = (b - t)/(b - a) 
and Wb{t) = (t - a) lib - a), and modified the approximation as 

N 

m « mm) ■.= f(f^!,)wait) + m)w,{t) + ^ r'^/wpsu, hxf^m, n) 

j=-N 

where 7"''' is defined by 7'"'[/](0 = f(t) - f(.f_^,^)Wa(t) - f{f^)wb{t). Throughout this paper, the 
formula © is called the generalized SE-Sinc approximation. 

2.3. (Generalized) DE-Sinc approximation and DE-Sinc indefinite integration 

Recently, instead of the SE transformation, the Double-Exponential (DE) transformation 

b — a ./^.. \ b + a 



b - a (K \ 

(/r°''(;c) = — — tanhl - smh(x)| + 



X = (p^'^t) = {i//°T\t) = arcsinh 



2 

— arctanh 

jt 



I2t-b-a\ 
\ b-a I 



has been employed by several authors [|3l |4l |8l [iS [l7|]. In this case, the formulas © and © are 
modified as 

N 



m « mm) := f{t^^^)wM + m)wi,{t) + J] T^^[mT)su,h)i(f>''\t)), (S) 

j=-N 

f(s)ds « Yj fi.tT){r^nih)j{iM{f\t)), (9) 

where tf = i/f°\jh) and r™[/](0 = f(t) - fit°%)Wa(t) - f{t''^)wb{t). Throughout this paper, the 
formulas ® and ® are called the generalized DE-Sinc approximation and the DE-Sinc indefinite 
integration, respectively. 

2.4. Convergence theorems 

Here let us introduce function spaces needed to state convergence theorems. 

Definition 1. Let ^ be a bounded and simply-connected domain (or Riemann surface). Then 
W°{S>) denotes the family of functions / analytic on & such that the norm ||/||H"(iS^) = sup^g^ |/(z)| 
is finite. 

Definition 2. Let a be a positive constant, and let ^ be a bounded and simply-connected domain 
(or Riemann surface) which satisfies (a, b) c ^. Then Lq,(^) denotes the family of functions 
/ e W°{&) for which there exists a constant L such that for all z in ^ 

\m\ < L\Q{z)\\ (10) 

where the function Q is defined by Q{z) = (z - a){b - z). 
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Definition 3. Let a be a constant with < a < 1, and let ^ be a domain with the same conditions 
as in Definition[2j Then Mq.(^) denotes the family of functions / 6 W°{^) for which there exists 
a constant M such that for all z in ^ 

\f{z)-f(a)\<M\z-ar, 
\f{b)-f{z)\<M\b-zV. 

In this paper, Qi is either of the following two domains: 

<A^=(^,) = {z = <A"(^) : ^ 6 ^rf} or r\%i) = {z = r\0 ■■ ^ e ^d), 

where ^j is the strip domain defined by ^^ = {^ 6 C : | Im ^| < J} for a positive constant d 
(see also Tanaka et al. [16, Figures 1 and 5] for the concrete shape of the domains). Convergence 
theorems for the generalized SE/DE-Sinc approximations are described as follows. 

Theorem 1 (Okayama |6, Theorem 3], see also Stenger [13]). Let f e M„(i/r^'^(^^)) /or d with 
< d < Ji, let N be a positive integer, and let h be selected by the formula 



, lid 



Then there exists a constant C which is independent of N, such that 

max \f(t) - r'^[fm\ <C^fNexp (- y/KdaN) . 

Theorem 2 (Okayama |6, Theorem 6]). Let f 6 Mai(l/°\%)) for d with <d < Jt/2, let N be 
a positive integer, and let h be selected by the formula 

log(2<«M) 

Then there exists a constant C which is independent ofN, such that 

-ndN 



max|/(0-n'[/](0|<Cexp 



a<t<b 1^ ^ ^ /v L^ . V . I - r y \og(2dN/a) ^ 
Convergence theorems for the SE/DE-Sinc indefinite integration have also been given as below. 

Theorem 3 (Okayama et al. iW, Theorem 2.9]). Let (fQ) e U,iil/'\^^)) for d with < d < n, 
let N be a positive integer, and let h be selected by the formula (fTTT) . Then there exists a constant 
C which is independent ofN, such that 

N 



max 

a<t<b 



< 



C exp (- ^IndaN^ . 



f{s)ds- Y, f(f;M''njh)j(j,h)(f\t)) 

j=-N 

Theorem 4 (Okayama et al. Q Theorem 2.16]). Let (fQ) e K(if/''^(^^)) for d with < d < 
jt/2, let N be a positive integer, and let h be selected by the formula (fT2l) . Then there exists a 
constant C which is independent ofN, such that 



max 

a<t<b 



N 



f{s)&s- Yj f{tf){r^njh)j{j,h){f\t)) 



j=-N 



AogCldn/a) ( -ndN 
< C exp ■ 



A^ ^ [\og(2dN/a) 



3. Numerical methods 

In this section, four numerical methods to be considered in this paper are described. First three 
methods are existing ones: the SE-Sinc-Nystrom method (Section 13.11) . the DE-Sinc-Nystrom 
method (Section [X2l) . and the SE-Sinc-coUocation method (Section [331) . Fourth one is the newly- 



proposed method: the DE-Sinc-coUocation method (Section [X4I) . 

3.1. SE-Sinc-Nystrom method 

As for the functions in Eq. ©, let y,(0 and gi{t) be each element of the vectors y{t) and g{t), 
respectively, and let kij{t) be (/, j)-th element of the matrix K{t). Assume the following conditions: 

(SEl) kijQ 6 Uir^^d)) {i=\,...,n, j=l,...,n), 

(SE2) giQ&l.Ar\%)){i=h...,n), 

(SE3) 3;;6H~(^^^(^rf))(/=l, ...,n), 

and define h as Eq. (fTTI) . Under those assumptions, the integral in Eq. ^ can be approximated by 
the SE-Sinc indefinite integration (|6]), and we have the new approximated equation: 

N 

i^^Xt) = r+Y, U(^f ) + K{tf)y'''\fl)] {^'YUhVih hXrm (13) 

In order to determine the approximated solution y^'^\ we have to obtain the unknown coefficients 
on the right hand side in Eq. (IT31) . i.e., 

F= = lyT\tt',), ..., yfV^), yf\f',), • • • , yfV^), • • • , y'fV^)]", 

which is a column vector of order (2 A^-l- l)-n (notice that n is the number of the system of equations, 
and A'^ is the number appearing in 2)- To this end, let us discretize Eq. (IT3l ) at (2A'^ + 1) sampling 
points: t = t^^ {i = -N, . . . , N), and derive the system of linear equations. Let us introduce some 
notation here. Let cTk = 1/2 + Si(jtA:)/jt, and let l''~^^ be a (2A^ + 1) x (2A^ + 1) matrix defined by 

4"'^ = [c^,-;], i,j = -N,...,N. 

Let In and /„ be identity matrices of order {2N +1) and n, respectively. Let D^ and Kfr be 
(2A^ + 1) X {2N + 1) diagonal matrices defined by 

D'^ = drngUr^n-Nh), ..., {i^'riNh)], 

KfJ = dmg[kij(tt'j,),...,kj(t%')], 

and let [Kff] be an n x n block of the matrices ^ff . Furthermore, let R and G^^ be column vectors 
of order (2A^ + 1) • n defined by 

R = [ri, . . . , ri, r2, . . . , r2, . . . , r„] , 

Then the system of equations to be solved is written as 

(/„ ®In- In ® {hI^-'^D'^}[K^f])Y''' = In ® {hf-'^D'^}G'^ + R, (14) 

where "(8>" denotes the Kronecker product. By solving the system (fT?l) . the approximated solution 
j*^^^ is determined by Eq. (IT3l) . This procedure is the SE-Sinc-Nystrom method. 
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3.2. DE-Sinc-Ny Strom method 

The important difference from the previous method is the variable transformation; the SE 
transformation is replaced with the DE transformation here. Assume the following conditions: 

(DEI) hjQ 6 U{r\^d)) {i=\,...,n, j=l,...,n), 

(DE2) giQ € l.„{r\%)) {i=\,...,n), 

(DE3) yi 6 H~(^°^(^rf)) {i=l,...,n), 

and define h as Eq. (fT^ . Under those assumptions, the integral in Eq. © can be approximated by 
the DE-Sinc indefinite integration ©J and we have the new approximated equation: 

N 

/''\t) = r+j] {gitf) + K(tf)/''\tf)] {r'njh)j(j, hxnt)). ds) 

j=-N 

In order to determine the approximated solution y^^\ we have to obtain the unknown coefficients: 

which is a column vector of order {2N + 1) • n. To this end, let us discretize Eq. (ITSl) at {2N + 1) 
sampling points: t = ?™ (/ = -A^, . . . , A^^), and derive the system of linear equations. Let D™ and 
Kfl^ be (IN + 1) X (2A^ + 1) diagonal matrices defined by 

D^^ = dmgUil^n'i-Nh), ..., {r^YiNh)], 
Kf[ = dmg[kj(l^l),...,kij{t''^)l 

and let [Kff] be an n x n block of the matrices Kff-. Furthermore, let G°'^ be a column vector of 
order (2N + I) ■ n defined by 

G''^ = [g,{t^%), ..., g,(tn g2(f%l • • • , g2(t7l ..., gnit'^N')]^. 

Then the system of linear equations to be solved is written as 

(/„ ®h-In® {hl'-'^D''^][K^fW^ = h® {hI^-'^D^^}G^^ + R. (16) 

By solving the system (fT6l) . the approximated solution j^^^ is determined by Eq. (fT5l) . This proce- 
dure is the DE-Sinc-Nystrom method. 

3.3. SE-Sinc-collocation method 

Stenger [13] developed the following SE-Sinc-collocation method independently of Nurmuham- 
mad et al. [5] (actually more than 10 years before), but below we find that it is strongly related to 
the SE-Sinc-Nystrom method described in Section 137X1 Assume the following conditions: 

(SEl) kjQ e Uir^i^d)) {i=h...,n, j=l,...,n), 

(SE2) giQ 6 Ui4f'\%)) (/=!,..., n), 

(SE4) yieM,{r\^d))ii=h...,n), 

and define h as Eq. (fTTI) . Let Y be the solution of the system of linear equations ([14] ). and let us 
write it as 

Y = [yi-N, yi,-N+l, ■ • ■ J yi,N, JI-N, JI-N+I^ • • • > y2,N^ • • • > yn,N\ • (17) 



Then the approximated solution j^^^CO = [yiit), . . ., yn(t)]^ is given by 

N 

yf(t) = yi,-NWa{t) + yuNWbit) + J] [yu - yi,-NWa(tf) - yi,NWb(tf)}S (j, h){f\t)), (18) 

j=-N 

for i = 1, . . . , n. This procedure is the SE-Sinc-coUocation method. 

3.4. DE-Sinc-collocation method (newly proposed) 

In view of Sections [3TT] and [X2l it is quite natural to replace the SE transformation with the DE 
transformation in the previous method. Assume the following conditions: 

(DEI) hjQ 6 U{r\^d)) (/ = 1, . . . , n, J = 1, . . . , n), 

(DE2) giQ 6 U{r\^d)) (^ = 1, . . . , n\ 

(DE4) y,eM,{r\^d)){i=h...,n\ 
and define h as Eq. (fT2l) . Let F be the solution of the system of linear equations (IT6] ). and let us 
write it as (IT7]) . Then the approximated solution j'^^^(0 = \yi{t), . . ., y„(0]^ is given by 

N 

yf\t) = y^,.NWait) + yuNWbit) + Yj iyiJ - yi,-NWa{tf) - yuNWb{tf)]S U, h)(cp''%t)), (19) 

j=-N 

for / = 1, . . . , n. This procedure is the DE-Sinc-collocation method. 

Remark 1. The assumptions on the solution j, i.e., (SE3), (DE3), (SE4), (DE4) seem to be hard to 
check, because y is an unknown function to be determined. In reality, however, those assumptions 
are unnecessary, because both (SE3) and (SE4) can be shown from the conditions (SEl) and (SE2), 
and both (DE3) and (DE4) can be shown from the conditions (DEI) and (DE2). To prove the facts 
is one of the main contributions of this paper, which is explained next (Theorems |6] and |7]). 

4. Theoretical results 

In this section. Theoretical results for the four methods in Section [3] are explained. The proofs 
are given in Section |6l 

4.1. Results on the regularity of the solution 

As described in Remark \T\ the condition on the solution y is assumed in each scheme. If the 
given problem ([T]) is a 'scalar' equation (n = 1), the following result has been known. 

Theorem 5 (Stenger et al. lll4 Theorem 2.3]). Let n = 1, and let the assumptions (SEl) and 
(SE2) be fulfilled. Then the initial-value problem ([T]) has a unique solution yi 6 'M^ai'^^i.^a))- 

This theorem shows the condition (SE4), and since M„(i^^''(^j)) c W°{il/\^J), the condition 
(SE3) is also shown. In this paper, the same result is shown in the case of a system of equations 
(for both SE and DE). 

Theorem 6. Let the assumptions (SEl) and (SE2) be fulfilled. Then the initial-value problem ([l) 
has a unique solution y with yt e Ma(if/^'^(^d)) for i = 1, . . . , n. 

Theorem 7. Let the assumptions (DEI) and (DE2) be fulfilled. Then the initial-value problem ^ 
has a unique solution y with yt e Ma(if/'^'^{^d)) for i = 1, . . . , n. 
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4.2. Results on convergence of the numerical solutions 

In the case of a 'scalar' equation, the convergence of the SE-Sinc-coUocation method is ana- 
lyzed as follows. In what follows, C denotes a constant independent of A^. 

Theorem 8 (Stenger UJ, pp. 446-447]). Let n = I, and let the assumptions (SEl) and (SE2) be 
fulfilled. Then, for all N sufficiently large, the system (fT4l) is uniquely solvable, and the error of 
the numerical solution yi ofEq. dTSl) is estimated as 

max \yxit) - yi{t)\ <CyfNexp (- ^IndaN) . 

This paper extends the result to a system of equations, and to the DE-Sinc-coUocation method. 

Theorem 9. Let the assumptions (SEl) and (SE2) be fulfilled. Then, for all N sufficiently large, 
the system (fT4l) is uniquely solvable, and the error of the numerical solution y'-^-' of Eq. (fT8l) is 
estimated as 

max (max \yi(t) - yf\t)\] < C ViVexp (- y/jtdaN) . 

i=l,...,n (a<t<b } ^ ' 

Theorem 10. Let the assumptions (DEI) and (DE2) be fulfilled. Then, for all N sufficiently large, 
the system (fT6l) is uniquely solvable, and the error of the numerical solution y'-^^ of Eq. (fT9l) is 
estimated as 

max (max \yi(t) - f/^m] < C exp ( "^ 

Furthermore, this paper also shows the convergence of the SE/DE-Sinc-Nystrom methods. 

Theorem 11. Let the assumptions (SEl) and (SE2) be fulfilled. Then, for all N sufficiently large, 
the system (fT4l) is uniquely solvable, and the error of the numerical solution y'^^'' of Eq. (fT3]) is 
estimated as 

max (max |y,(0 - 3;;"^)! 1 < C exp (- y/jidaN) . 

i=l,...,n {a<t<b I ^ ' 

Theorem 12. Let the assumptions (DEI) and (DE2) be fulfilled. Then, for all N sufficiently large, 
the system (fT6l) is uniquely solvable, and the error of the numerical solution y'^^'' of Eq. (flSl) is 
estimated as 

log(2dn/a) ( -ndN 



max max ^t) - yf\t)\ < C ^\' ' exp 

i=\,...,n \a<t<b ) N 



log(2dN/a) 



4.3. Discussion about the performance 

In view of the convergence rates shown above, the DE-Sinc-Nystrom method seems to be 
the best, and this was then followed by the DE-Sinc-coUocation method, the SE-Sinc-Nystrom 
method, and the SE-Sinc-coUocation method. However, the DE-Sinc-coUocation method (the sec- 
ond one) can be considered as the best, or at least as useful as the DE-Sinc-Nystrom method, 
for the following reasons. Firstly, the difference of convergence between the DE-Sinc-Nystrom 
method and the DE-Sinc-coUocation method is quite small, and actually it is almost indistinguish- 
able in the numerical experiments (see Figures [l3-|2]in Section [5]). Secondly, compared to the 
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the approximate solution of the DE-Sinc-collocation method y'-^-' (Eq. (fT9l)). that of the DE-Sinc- 
Nystrom method y^^^ (Eq. (flSl) ) has time-consuming terms to evaluate. All of the basis functions 
in y'-^^ are elementary functions, whereas the basis functions J(j, h) in y^'^* includes the special 
function Si(x). Furthermore, j^^^ can be computed with 0{nN), but y^^^ needs 0{n^N) because a 
matrix- vector product is included in y^'^K Therefore, from the viewpoint of the computational cost, 
the DE-Sinc-collocation method is better than the DE-Sinc-Nystrom method (see also Table [U). 

5. Numerical results 

In this section, numerical examples of the SE/DE-Sinc-Nystrom methods and the SE/DE-Sinc- 
coUocation methods are presented. The computation was done on Mac OS X 10.6, Mac Pro two 
2.93 GHz 6-Core Intel Xeon with 32 GB DDR3 ECC SDRAM. The computation programs were 
implemented in C-I--I- with double-precision floating-point arithmetic, and compiled by GCC 4.0.1 
with no optimization. The linear systems (fT4l) and (fT6l) are solved by using the LU decomposition. 
In what follows, jt_ denotes an arbitrary positive number less than Jt, and it was set as Jt_ = 3.14 
in actual computation. Firstly, let us consider the following two examples. 



Example 1. Consider the following initial value problem (the Halm equation [Il2ll ) over the inter- 
val [0, 1]: 

il+x^)Yit)-2y = 0, yiO) = 0, /(O) = 1, 

which is equivalent to the system 

y[(t) = yiit), ym = 0, 

y^^^^ = n 1 2^2 >^l(^)' 3^2(0) = 1, 

(1-1- x^Y 

whose solution is yi{t) = y/l + x^ sinh(arctan x), y2(t) = y[{t). 
Example 2. Consider the following initial value problem over the interval [0, 2]: 

y[(t) = -y,it) + -^y2(t), ym = 0, 

y'lit) = — ?yi(t), yiiO) = 1, 

whose solution is yi(t) = V^exp(-0, yiit) = exp{-t). 

As for ExamplelU the conditions (SEl) and (SE2) are satisfied with a = 1 and d = 3jt_/4. In 
the DE case, let us set p = Jt_/(2 log 2) and 



q= ^{l + Ip^ + V(l + 7p2)2 + i6py]j2, 



and furthermore set x = - (1 - q) /(4-p), y = 3(1 - (l/q)) /4, and J_ = aicsiniy/ yjx^ + y^). Then, 
the conditions (DEI) and (DE2) are satisfied with a = 1 and d = d^. As for Example |2l which is a 
harder example because of the singularity at the origin, (SEl) and (SE2) are satisfied with a = 1/2 
and d = Jt„, and (DEI) and (DE2) are satisfied with a = 1/2 and d = jt„/2. The numerical errors 
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Figure 1 : Errors in Example [T] 
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Table 1: Computation times and A^ needed to obtain 10 '^ accuracy in Example |2] 
SE-Sinc-Nystrom SE-Sinc-collocation DE-Sinc-Nystrom DE-Sinc-collocation 



87 
0.281 



87 
0.137 



31 
0.107 



31 
0.050 



are plotted in Figures [Hand |2l respectively. In the graphs, "maximum error" denotes the maximum 
absolute error at 999 equally-spaced points (say f/) on the interval [a, b], i.e., 



maximum error = max j max |y,(?/) - y,(?/ 

r=l,...,nW=l,...,999 



where yj means each numerical solution. From both figures, we can confirm the results of The- 
orems l9l-fT2l More precisely as for the (newly -proposed) DE-Sinc-collocation method, its con- 
vergence rate is actually much higher than that of the SE-Sinc-collocation method. As described 
in Section 14.31 although the theoretical rate of the DE-Sinc-Nystrom method is a bit higher than 
that of the DE-Sinc-collocation method, both rates are almost indistinguishable in the numerical 
results. Moreover, as seen in Table [1] the DE-Sinc-Nystrom methods needs times twice as much 
as the DE-Sinc-collocation method to obtain 10"^ accuracy (the same applies in the SE case). At 
least from the result, we can conclude that the DE-Sinc-collocation method is the most efficient. 

In the examples above, all the assumptions (SEl), (SE2), (DEI), and (DE2) are satisfied with 
some a and d. Let us have a look at another case here. 



Example 3. Set a function F as F{t) = Vcos(4 arctanh t) + cosh(jt), and consider the following 
initial value problem over the interval [-1, 1]: 



y\(t) = 



2[tF^(t) + sin(4 arctanh t)] 

m 



yiit). 



2[tF^{t) + sin(4 arctanh t)] 



Ji(-1) = 0, 

yii-i) = 1, 



whose solution is yi{t) = sin[(l - f)F(t)], j2(0 = cos[(l - f)Fit)]. 
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Figure 3: Errors in Example |3] 



This is a quite hard example to solve numerically, due to the bad behavior of F at f = ±1 
(non-regular points are densely distributed around the endpoints). Fortunately, the assumptions 
(SEl) and (SE2) are satisfied with a = I and d = Jt_/2, but (DEI) and (DE2) are not satisfied with 
any d > (we easily see a = 1, though). Therefore, Theorems [TO] and [T2] cannot be used in this 
case. However, according to the recent result [lllfl . even in such a case, DE's methods may achieve 
the same convergence rate with that of SE, by setting d = arcsin(JsE/^), where Jse denotes SE's d. 
We can in fact observe it in Figure |3j DE's methods seem to converge with the similar rate to that 
of SE. Since the computational cost is the same as that of the previous examples, we can consider 
that the DE-Sinc collocation method still keeps the lead even in this case. 



6. Proofs 

6.1. Proofs on the regularity of the solution 

The idea here is to apply the standard contraction mapping theorem, which holds not only in 
the scalar case but also in the case of a system of equations. SetX = {H°°(^)}'' and Y = {Mq,(^)}", 
and define ||/||x = niax,=i ,„{||yi||H~(st)}. The goal is to show j e Y, but it is not easy because Y is 
not a Banach space. For this reason, firstly j e X is shown (X is a Banach space), and by using the 

f{s) ds, 
andT:X^Xas 



^[/](0 



•Ja 



K(s)f{s)ds, 



where K satisfies the assumption (SEl) or (DEI). If the operator is multiplied repeatedly, it be- 
comes a contraction map. 

Lemma 13. Let the assumption (SEl) be fulfilled. Then it holds for all positive integers m and 

z 6 (/^""(^rf) that 



l^'"[/](z)l < 



{nL(b - ay"-'ci B(i^i(x), a, a)}" 



ml 



■Il/llx[l, 1, ...,lf, 
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where x = Ke\il/^{z)\, if/i(x) = (tanh(x/2) + l)/2, B{x,a,/3) is the incomplete beta function, L is 
the constant in Eq. (fTOl) . and C\ is a constant depending only on d. 

Lemma 14. Let the assumption (DEI) be fulfilled. Then it holds for all positive integers m and 
z e iA™(^d) that 

m! 
where x = Re[i/r°'^(z)], iA2(-'c) = (tanh(7t sinh(x)/2) + l)/2, L is the constant in Eq. (fTOl) . anJ C2 /^ a 
constant depending only on d. 

These lemmas are straightforward extension from the existing ones [i9|, Lemmas 5.4 and 5.6], and 
the proofs are omitted. Then in both cases it holds that 

m! 
and thus for sufficiently large m, 'V'" is a contraction map, from which we have the next theorem. 

Theorem 15. Let the assumptions (SEl) and (SE2) be fulfilled. Then Eq. ([2]) has a unique solution 
yeX, i.e., yt e H~(^^^(^^))/or i=l,...,n. 

Theorem 16. Let the assumptions (DEI) and (DE2) be fulfilled. Then Eq. ([2]) has a unique solu- 
tion yeX, i.e., yi e H~(i/^™(^rf))/or i=l, ...,n. 

Proof. Before applying the contraction mapping theorem, the only thing we have to show is 
(SE2)/(DE2) ^JgeX, which is done by Lemmas [H] and [13 (since M„(^) c H"'(^)). 

The next lemma is a result for SE, which completes the proof of Theorem [T5l 

Lemma 17 (Stenger iQ Theorem 4.1.3]). Let gQ e Lq.(0-"'(^j)), and set q(t) = f'g(s)ds. 
Then q e M„(i/r^^(^^)). 

In the case of DE (for Theorem [16]), we need the next lemma. 

Lemma 18 (Okayama et al. IlZI, Lemma A.4]). For x eR andy e (-tc/2, tc/2), it holds that 

1 , /Jtcosj . , \ 1 1 , /^ . , . . .\ 1 
tj^iix) := - tanh ^ — - — smh ^J + 2 - 2 \2 ^ *^^ ^^^7^2 " 

Using Lemma [TSl we show the following lemma (DE version of Lemma [T71). 
Lemma 19. Let gQ e l.^iip^'^^d)), and set q{t) = f^g{s)ds. Then q e Ma{if/'"^(%)). 
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Proof. Clearly q 6 W°{ifj°^{2id)) holds. Let us show the Holder continuity at the endpoint a 
(showing it at b is omitted, because it is quite similar). Putting n=l,/=l,^ = ^in Lemma [141 
we have 

<L{b-af"-'c2'R{ilf2{x),a,a). 






g{w) dw 



Then it holds that 



\q{z) - q{a)\ = 



Jg{w) dw 
a 

Furthermore, from Lemma [18] it holds that 

(b - a)i//2(x) < 



< Lib - af"-^CdB(il/2(x),a,a) < Lib - af"-^CdBil,a,a){il/2ix)f. 



b - a , (^ . , . \ b - a 

—— tanhl - smhix + iy)\ + — — 



= |i/r™(jc + iy)-a| = \z-a\. 
Thus there exists a constant L such that \qiz) - qia)\ < L\z - af for all z 6 il/°^i^d)- D 



Now showing j e X is finished. For y e Y, what is left is to show the Holder continuity. In 
view of the right hand side of Eq. (O, clearly y (on the left hand side) is Holder continuous of 
a-order. Hence Theorems [6] and [7] are established. 

6.2. Proofs on convergence of the numerical solutions 
6.2. L SE/DE-Sinc-Nystrom method 

Firstly, the SE-Sinc-Nystrom method is considered. Notice that in this subsection, set C = 
{Ci[a, ^])}", and all operators here are discussed on this function space. Let us introduce the 
operator ^^'^, which is an approximation of ^, as 



jNiim = 2 ntf){r'Yijh)Jij,h)irm 



-N 



and "V^ as 'y^[/](0 = SF^[Kf]it). Then consider the following three equations: 



iJ-<y)y = r + Jg 

iln ^In- In ® {hl'^-'^D'^mfW^^ = /„ {hI^-'^D'^}G'' + R 



(Eq. @), 

(Eq. m), 

(Eq. ([H)). 

Using the standard arguments (e.g., see [|9l Lemma 6.1]), we can see the unique solvability of 
Eq. ([14]) is equivalent to that of Eq. ([T3]) . If the unique solvability of Eq. ([T3]) is shown, i.e., 
(J - "V^)"' exists, we have 

= (J - ^'j^r' {(r + Jg + ^y) - ^'^y - r - J'^g] 

= il- ^^^)-i [iJg - J'^'g) + i<Vy - ^'^'y)] , 

and finally using Theorem [3] the desired error estimate (Theorem[lT]) is obtained. Therefore, what 
is left is to show the existence and boundedness of (J - "V^)"' . For the purpose, the next theorem 
is useful. 

Theorem 20 (Atkinson [1, Theorem 4.1.1]). Assume the following four conditions: 

14 



1 . Operators K and Xn are bounded operators on C to C. 

2. The operator {I - X) : C ^ C has a bounded inverse (J - X)~^ : C ^ C. 

3. The operator Xn is compact on C. 

4. The following inequality holds: 

\\{X - Xn)Xn\\uC,C) < 7777 TJ— 77] • 

||(i -A) MIX(C,C) 

Then {I - X„)~^ exists as a bounded operator on C to C, with 

M.j- wx-l|| 1 + ||(J - ^r'|lx(C.C)ll^«llx(C.C) 

"^^"'^"^ "^^^^^^-l-||(J-^)-|lx,c,c)ll(^-^.)^Jx(CC)- ^''^ 

We need to show the four conditions in this theorem as /Y = "V and X„ = "V^. The first 
condition clearly holds, and the second condition is known as a classical result. The third condition 
immediately follows from the Arzela-Ascoli theorem. The fourth condition is shown by the next 
lemma, which is straightforward extension from the existing one [^ Lemma 6.5]. 

Lemma 21. Let the assumption (SEl) be fulfilled. Then there exists a constant C independent of 
N such that 

Furthermore, ll'y^llx(c.c) is uniformly bounded, since "V^/ converges to "V/ for any / 6 C. 
Thus, from Eq. (|20] ). we obtain the desired result: (J - "V^)"' exists and uniformly bounded for 
all sufliciently large N. This completes the proof of Theorem [TTI (the SE-Sinc-Nystrom method). 

The proof for the DE-Sinc-Nystrom method goes on in exactly the same way. Let us introduce 
the operator ^™ as 

N 

jN^wm = Yj fitf){r"nih)j{fh){f\t)\ 

j=-N 

and "V™ as "V^L/ICO = J^N^[Kf](t). The diff'erence from the SE is the next lemma, which is also 
straightforward extension from the existing one [9, Lemma 6.9]. 

Lemma 22. Let the assumption (DEI) be fulfilled. Then there exists a constant C independent of 
N such that 

IK V - Vf^)VN llx(c,c) s c I — — 

This completes the proof of Theorem [T2] (the DE-Sinc-Nystrom method). 

6.2.2. SE/DE-Sinc-collocation method 

Let us consider the SE-Sinc-coUocation method first. Notice the relation y-^V) = ^^[j;^^](0> 
where y. is the solution of the SE-Sinc-Nystrom method (see Eq. (fT3T )). and y. is the solution 
of the SE-Sinc-coUocation method (see Eq. (fTSl) ). Then we have 

Wyi-y^Wc < Wyi-KyiWc + WKUcoWyi-yTWc (21) 
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The first term on the right hand side can be estimated by Theorem [TJ On the second term, use 
Theorem nn for ||y,- - jI'^^Hc, and use the next lemma to obtain ||P^||£(c,c) ^ C log{N +1). 



Lemma 23 (Stenger ||13|, p. 142]). Let h > 0. Then it holds that 

supV |50;/i)(^)|<-(3 + logA^). 

This completes the proof of Theorem |9] (the SE-Sinc-coUocation method). 

The proof for the DE-Sinc-coUocation method goes on in exactly the same way. By using the 
relation y^^\t) = 'P^[y''f^h(t), we have the similar inequality as Eq. (I2TI) (just replace SE with DE). 
By estimating each term via Theorems l2] and \T2\ and Lemma l23l Theorem [TOl is established. 
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